Connectivity-dependent properties of diluted sytems 
in a transfer-matrix description 
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We introduce a new approach to connectivity-dependent properties of diluted systems, which is 
based on the transfer-matrix formulation of the percolation problem. It simultaneously incorporates 
the connective properties reflected in non-zero matrix elements and allows one to use standard 
random-matrix multiplication techniques. Thus it is possible to investigate physical processes on the 
percolation structure with the high efficiency and precision characteristic of transfer-matrix methods, 
while avoiding disconnections. The method is illustrated for two-dimensional site percolation by 
calculating (i) the critical correlation length along the strip, and the finite-size longitudinal DC 
conductivity: (ii) at the percolation threshold, and (Hi) very near the pure-system limit. 
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The transfer-matrix (TM) approach to percolation was 
pioneered by Derrida and coworkers [Q . By analysing the 
possible combinations of adjacent column states made up 
of occupied and unoccupied sites (or bonds), and the al- 
lowed connections among the latter and to the arbitrary 
origin of a two-dimensional (d = 2) strip, it was possible 
to write the TM on the basis of such column states. The 
key element in this formulation was the fact that, from 
the very structure of the TM, repeated multiplication is 
tantamount to the simultaneous generation of all pos- 
sible connected configurations that span the strip, each 
with its proper probabilistic weight. Thus the probabil- 
ity of connection to the origin, whose exponential decay 
is governed by the correlation length, is asymptotically 
given exactly by the largest eigenvalue of the TM. The 
correlation length could then be used in a phcnomcno- 
logical renormalisation calculation g], which gave very 
accurate results for critical parameters such as the per- 
colation threshold and correlation-length exponent v. It 
was not clear, however, how one could take advantage of 
such a direct and elegant scheme to investigate properties 
other than the decay of the probability of connection to 
the origin. The obvious alternative, of building up suc- 
cessive columns by occupying (or not) each individual 
site independently, runs into the problem of disconnec- 
tions, which is severely aggravated on a strip geometry. 
Up to now, the usual solution has been to generate con- 
figurations site by site, and study quantities which do not 
depend on keeping connectivity along the strip, e.g. the 
moments of the distribution of clusters H for percola- 
tion in d = 2 and 3. A clever way to get round the effects 
of disconnections for random resistor-insulator networks 
at percolation was introduced Q by generating individ- 
ual elements on long strips (or bars, in d = 3) with free 
edges. By imposing a fixed voltage drop across the strip, 
it was possible to invoke TM concepts in a step-by-step 



evaluation of the transverse conductivity, for which lon- 
gitudinal disconnections of the resistor structure are ir- 
relevant. Finally, in superconductor- resistor networks at 
the percolation threshold of superconducting elements, 
disconnections are in fact responsible for the quantity of 
interest, which is the residual finite-size resistivity; one 
can then establish periodic boundary conditions across 
the strip (in order to minimise finite-width effects) and 
estimate the longitudinal resistivity 

Here we introduce a scheme which preserves the con- 
nected structure of the percolating cluster as one sweeps 
along the strip, and at the same time relies on standard 
ideas of random-matrix multiplication. The latter fea- 
ture implies that any physical quantity, in addition to 
connection probability, can be sampled along the strip 
through insertion of its corresponding local realisation. 
This opens the way e.g. to the straightforward treat- 
ment of spin-spin correlations in dilute magnets 0, for 
which so far only approximate TM treatments, relying 
on plausible but essentially uncontrollable assumptions, 
have been available || . In the new scheme we enumerate 
the set of all allowed column combinations, according to 
the original TM procedure jjj, and then build the strip 
one full column at a time, by picking a given column's 
successor at random but only among those columns that 
are allowed by the connectivity rules (i.e., which have a 
non-zero TM element linking them to the immediate pre- 
decessor). With the proper assignment of probabilistic 
weights, as explained below, this procedure is equivalent 
to the sampling of connected configurations implicit in 
the iteration of the TM. 

In what follows, we first expose the basic concepts of 
the method; then the decay of correlations is calculated 
and shown to reproduce the results given by diagonal- 
isation of the TM. Next we apply the method to the 
longitudinal conductivity of a diluted resistor-insulator 
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network. For this particular process high-accuracy re- 
sults exist, together with some exact ones, which provide 
a test of the method. At the percolation point, the new 
method produces estimates of the conductivity exponent 
which compare very favourably with those existing in the 
literature. Near the pure- network limit, we obtain the 
corrections to the conductivity to first order in defect 
concentration, in excellent numerical agreement with an- 
alytical results. Finally we point to possible extensions 
and generalisations of the present approach. 

We consider a strip of a square lattice of width L sites, 
with periodic boundary conditions across, on which sites 
may be independently occupied (unoccupied) with prob- 
ability p (1 — p). As explained in Ref. [Q, one builds all 
possible column configurations in which at least one site 
is occupied and linked to the origin (assumed to be at the 
far left, say); other sites may be either occupied and con- 
nected, occupied and disconnected, or unoccupied. The 
TM element Tij between column configurations i and j 
is non-zero only if column j is allowed to succeed column 
i (meaning: connection to the origin must be preserved, 
and illogical situations, such as an occupied and discon- 
nected site being preceded by a connected one on the 
same row, must not occur). One has 



P N ( P ) ~exp(-iV/a(p)) 



(4) 



T ij =p N '{\-p) 



L — Nj 



(1) 



where Nj is the number of occupied sites in j. Our pro- 
cedure then goes as follows. Assume that the strip has 
been built up to a column whose configuration is i. Call 
= ■ ■ -jMi the set of all j's (a total of Mi) 

allowed to succeed a given i. A segment of length 



U = 



Lij — 



p Nj (l-p) 



{3(0} 



L — Nj 



(2) 



represents the total (conditional) probability to have a 
connected configuration succeeding column i. Drawing a 
random number < e < 1 from a uniform distribution, 
the next column configuration is chosen to be ji such 
that 



3=31 



Li a < eLi < 



3i 

^2 Li 3 
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This ensures that the allowed connected configurations 
come up with their proper corresponding probabilities. 
One can then proceed and generate a column to succeed 
jo, and so on iteratively. The only information to be kept 
in store throughout the process is the same as that used 
in the standard TM formalism: the location and indices 
(column occupancy numbers) of non-zero TM elements. 

We now show that for a strip of width L, the scheme 
described above gives the same correlation length, 
as that obtained from diagonalisation of the TM. One 
defines £l (p) through the exponential decay of the prob- 
ability of connection between columns and A, Pjv(p): 



As the process described here is a sequential sampling 
one (as opposed to the parallel updating of mutually ex- 
cluding paths in configuration space, which takes place 
in the iteration of the TM), one must consider the appro- 
priate quantities to analyse. At each step, when column 
configuration jj is chosen among Mi to succeed i, one 
is probing one branch of a tree structure in the space 
of column configurations, and discarding Mi — 1 others. 
In order to deal with this, and produce unbiased sam- 
ples, the standard procedure is the weighting of steps 
introduced in early simulations of self-avoiding walks M . 
When considering the allowed moves from a point i to 
the next (here: in configuration space, as opposed to real 
space in Ref. [|| ) one generates a weight Wi proportional 
to the total probability of moving out from that point. 
In the present case, Wi = Li of Eq. (||). It is easy to 
see that W% is properly normalised, as the denominator 
is the sum of probabilities of all possible succeeeding con- 
figurations, not only connected ones, and that is unity. 
The total weight of a given N — 1- step walk (spanning 
A points) is the product W1W2 ■ ■ ■ Wn, to be denoted 
by Wat. The quantity whose variation with distance is 
to be studied, in the present context of sequential sam- 
pling, is WV; on universality grounds, it is expected to 
decay with the same correlation length as Pn of Eq. (Q) . 

For strips of widths L = 4, ...11 at p = 0.592745, 
the best numerical estimate of the appropriate perco- 
lation threshold to our knowledge |10|, we have gener- 
ated large numbers (N samp ) of independent connected 
configurations between the origin and column Nq. The 
weights of the N samp configurations up to respectively 
columns Ao — A and Nq were summed to produce the es- 
timates W No - A = {Ek=T F W ( No-A)/ N samP and Wjv„ = 
(J2k=i" P )/iV som p . An estimate of the correlation 
length is then given by 



A _ ln pV^_ 
£ \Wn -a 



(•5) 



Finally, we have repeated the process n s times with dis- 
tinct random- number sequences, in order to estimate 
fluctuations. We used N samp = 10 5 , A = 20, A = 10, 
n s = 100. Our results are displayed in Table 1, in the 
form of estimates for the critical decay-of-correlations ex- 
ponent r), through the identity 77 = L/n^(p c ) given by 
conformal invariance |ll|l . These are to be compared with 
those, also in Table 1, obtained from the largest eigen- 
value of the TM. The values of Ao and A were chosen 
bearing in mind that, both from general finite-size scal- 
ing ideas and from previous results for percolation [JIJ , it 
is known that the correlation length at criticality must 
be of order L. The convergence of finite- width results 
towards the value given by conformal invariance [ pT| , 
1] = 5/24 = 0.208333... has been investigated else- 
where |l2| ]. For our present purposes the relevant com- 
parison is between the columns of Table 1, which shows 
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the soundness of the proposed scheme. As expected from 
the theory of normal distributions, fluctuations shrink 
with (N samp ) , because it is there that the accumu- 
lation of sample weights leading to self-averaging takes 
place, but remain approximately constant with n s . The 
amount of computational time involved is linearly pro- 
portional to n s x N samp , thus one can produce more ac- 
curate results by increasing N samp while reducing n s to 
one; in this limit, the width of the error bar for the single 
(presumably very precise) central estimate can be extrap- 
olated from those obtained for large n s and correspond- 
ingly smaller N samp , via the (N samp ) dependence. 
Since our goal here is to demonstrate the feasibility of 
the proposed approach, rather than refining numerical 
values, we did not pursue this line systematically. 

We now show results for finite-size conductivity o\l(Pc) 
at the percolation threshold of a resistor-insulator net- 
work. From finite-size scaling, this is expected to vary 
with strip width as 

vl(Pc) ~ L-V" , (6) 

where the best available estimate for the exponent is 
t/v = 0.9745 ± 0.0015 [§]. We have generated samples 
of site-diluted resistor networks (where a bond is a re- 
sistor if it connects two occupied sites, and an insulator 
otherwise), according to the procedure delineated above. 
Now, since the quantity to sample (average conductiv- 
ity per bond) is naturally accumulated as one proceeds 
along the strip (instead of decaying exponentially, as is 
the case with connection probabilities) , one does not have 
to be concerned with weights: it suffices to generate very 
long samples, and each column configuration will come 
up with its good weight. 

Conductivities have been calculated by Fogelholm's 
node deletion algorithm |l3|. This is very efficient on 
a strip geometry since it depends on keeping track only 
of at most L(L — l)/2 links among sites, plus 2L links to 
the origin. 

Table 2 shows our data for L — 3-11 where for each 
strip width 100 independent samples, each of length 10 5 
columns, were generated. Error bars reflect deviations 
among different samples. A least-squares fit to a log-log 
plot of the data in Table 2 gives t/v = 1.005±0.002, with 
an accumulated \ 2 P er degree of freedom (DOF) = 1.0, 
an estimate which is 3% above the accepted value 
with apparently non-overlapping error bars. Before ac- 
cepting this at face value, some remarks are in order. 

Firstly, we perform a similar fit to the resistivity (pl) 
data on the site superconductivity problem, shown in Ta- 
ble I of Ref . U . In d = 2 the superconducting exponent 
s/v in pl ~ L~ s l v is the same as t/v, by duality (see 
Refs. f|||). Also, the same (periodic) boundary condi- 
tions across the strip were used there, as opposed to free 
ones for previous conductivity studies Q . This is impor- 
tant when one wishes to compare purely finite-size effects 
between two sets of data. Those authors simulated strips 
10 4 times as long as we did, thus the fact that their error 



bars are two orders of magnitude smaller than our own 
indicates that both methods have the same intrinsic ac- 
curacy. We then turn to comparison of systematic errors. 
For the same range 3 < L < 11, the data of Ref. p] give 
t/v ~ 0.913, though with a very large x 2 /DOF which 
partly reflects the greater accuracy of individual data in 
Ref. |(|, as well as the need to take corrections to scal- 
ing into account. Assuming a power-law correction with 
exponent w, always for 3 < L < 11, fits of plL 1 ^ to 
a + b/L u for 0.94 < t/v < 1.02 show that x 2 / D F in- 
deed has a minimum value at t/v ~ 0.982-0.987 when 
u> is kept constant at 1.2 - 1.4 (see Ref. Q). The am- 
plitude b varies monotonically between —0.1 and —0.6. 
A similar analysis of our own data shows that x 2 /DOF 
has a gentle maximum at t/v ~ 0.95 and a minimum 
at t/v ~ 1.025. The amplitude b starts from 0.44 at 
t/v = 0.94 and decreases monotonically, crossing zero at 
t/v ~ 1.005. Varying w along a wider interval, between 
1 and 2.5, does not produce any significant change. 

Thus, for similar strip widths and comparable amounts 
of computational effort, our method generates data of 
comparable quality to other authors'. It seems that for 
critical conductivity studies in two dimensions one has 
to reach very large widths, of order 40 sites Q-^j], be- 
fore asymptotic behaviour sets in. While this could be 
done in Refs. the nature of the present algorithm 

is such that the exponential growth, with strip width, 
of the number of configurations to be stored is the main 
obstacle to going further than L — 11. However, from 
past experience |8[ we expect such upper limit not to be 
as stringent for e.g. diluted magnets. 

We have also studied resistor networks for very low 
impurity (insulator) concentrations (1 — p), where it has 
been predicted fl4|| that conductivity must vary as 

a(p)/a(l) = l-7r(l-p)+7r(l-p) 2 /2 . (7) 

By using finite-size considerations pertaining to low con- 
centrations of impurities in a cylindrical geometry, we 
have shown that the finite-size conductivity is 

a L (p)=a 00 (p) + a(p)/L 2 + 0(l/L 4 ) , (8) 

where (Too(p) is given by Eq. ([?]) and a(p) = — (1— p)n 3 /6. 
Our data for the normalised conductivity at p — 0.999 
(where Eq. (p|) gives CToo(p) = 0.99686... and a(p) — 
-5.17 x 1Q-3) and 3 < L < 11 are shown in Table 3, 
where for each strip width 100 independent samples, each 
of length 10 5 columns, were generated. Error bars reflect 
deviations among different samples. A least-squares fit of 
our finite-size data gives Coo = 0.99685 ± 0.00005, where 
the small error bar reflects the fit's overall smoothness, 
and a = (—6. Oil. 7) x 10 -3 , in very good agreement with 
the theoretical prediction. 

We have proposed and illustrated a straightforward 
scheme for diluted systems, in which a transfer-matrix 
approach can be implemented without giving rise to lon- 
gitudinal disconnections along a strip. Previous treat- 
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merits either were restricted to the calculation of the de- 
cay of connection probability |jj , or could be carried out 
independent of disconnections, owing to particular geo- 
metric [|] or physical features, or else were forced 
to rely on essentially uncontrollable assumptions on the 
commutation of TMs associated to distinct dilution con- 
figurations @] . Extensions of the present work to dilute 
magnets are now being considered. Further appli- 
cations would be to the anomalous thermal behaviour of 
Fe(llO) submonolayers on W(110) |lq| , and to frustrated 
percolation |L6]], a problem related to glass- formation 
processes. 
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TABLE I. Estimates of n = L/n£(p c ) at p = 0.592745; 
averages over n s = 100 distinct sequences, each of 



N sar . 



10 accumulated weights . 



L 

4 
5 
6 
7 
8 
9 

10 
11 



77 (this work) 

0.21255 ±0.00020 
0.21142 ±0.00026 
0.21069 ±0.00029 
0.21016 ±0.00039 
0.20979 ± 0.00043 
0.20954 ±0.00044 
0.20930 ±0.00047 
0.20915 ±0.00055 



X] (TM) 

0.2125576128 
0.2114673276 
0.2107370714 
0.2102232886 
0.2098564767 
0.2095868033 
0.2093833099 
0.2092261631 
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TABLE II. Estimates of o L {Pc) at p = 0.592745; averages 
over n s = 100 distinct strips, each of length 10 5 columns . 



L 




3 


0.34435 ± 0.00063 


4 


0.25931 ±0.00054 


5 


0.20718 ±0.00053 


6 


0.17241 ±0.00047 


7 


0.14749 ±0.00045 


8 


0.12881 ±0.00038 


9 


0.11431 ±0.00042 


10 


0.10271 ±0.00035 


11 


0.09323 ±0.00033 



TABLE III. Estimates of a L (p) at p — 0.999; averages 
over n s — 100 distinct strips, each of length 10 5 columns. 
Extr.: extrapolation against 1/L 2 (see text). Expected: 
Eq. ®. 



L a L 

3 0.99618 ±0.00019 

4 0.99648 ±0.00016 

5 0.99661 ±0.00013 

6 0.99669 ±0.00012 

7 0.99673 ±0.00011 

8 0.996753 ±0.000099 

9 0.996775 ±0.000099 

10 0.996791 ±0.000095 

11 0.996804 ±0.000092 
Extr. 0.99685 ±0.00005 

Expected 0.99686 . . . 
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